clear all; close all; clc;
location = strcat(pwd,'\');
varout = strcat(location,'data\');
path(path,strcat(location,'function_calls'))

% Parameters
A = 22; % leverage
K = 0.05; % fixed cost
ann_alpha = 3:0.3:9; % true alpha [3 to 9]
ann_sigma = 11:0.8:19; % signal precision [aum+returns]
gamma = 12; % cost convexity parameter
rho = 0.05/12; % intertemporal discount rate

maxit= 300; % for hbj iteration
crit = 10^(-6);

% State space grid
I=101; % p bounded [0,1]
J=101; % q bounded [0,1]

%% Vary alpha in this window

for ali=1:size(ann_alpha,2)
    
    alpha = ann_alpha(ali)/(12*100);
    
    ind = (size(ann_sigma,2)+1)/2;
    sigma = ann_sigma(ind)/(sqrt(12)*100); % mid sigma
    
    [v,inten,pi,p,q,sigma,alpha] = HJB_diffusion_implicit(A,K,alpha,...
        sigma,gamma,rho,maxit,crit,I,J);
       
    % Save value function and other details (given i,j)
    % note strcat can't handle decimals, multiple parameter values by 10
    sub = strcat('_',num2str(ann_alpha(ali)*10),'_',num2str(ann_sigma(ind)*10)); % renames
    
    S.(strcat('v',sub)) = v;
    S.(strcat('inten',sub)) = inten;
    S.(strcat('pi',sub)) = pi;
    S.('p') = p;
    S.('q') = q;
    S.(strcat('var',sub)) = sigma;
    S.(strcat('alpha',sub)) = alpha;

    save(strcat(varout,'vf','_',num2str(ann_alpha(ali)*10),'_',num2str(ann_sigma(ind)*10),'.mat'),'-struct','S');    
end

%% Vary sigma in this window

for sigj=1:size(ann_sigma,2)
    
    ind = (size(ann_alpha,2)+1)/2;
    alpha = ann_alpha(ind)/(12*100);
    sigma = ann_sigma(sigj)/(sqrt(12)*100); % mid sigma

    [v,inten,pi,p,q,sigma,alpha] = HJB_diffusion_implicit(A,K,alpha,...
        sigma,gamma,rho,maxit,crit,I,J);
       
    % Save value function and other details (given i,j)
    % note strcat can't handle decimals, multiple parameter values by 10
    sub = strcat('_',num2str(ann_alpha(ind)*10),'_',num2str(ann_sigma(sigj)*10)); % renames
    
    S.(strcat('v',sub)) = v;
    S.(strcat('inten',sub)) = inten;
    S.(strcat('pi',sub)) = pi;
    S.('p') = p;
    S.('q') = q;
    S.(strcat('var',sub)) = sigma;
    S.(strcat('alpha',sub)) = alpha;

    save(strcat(varout,'vf','_',num2str(ann_alpha(ind)*10),'_',num2str(ann_sigma(sigj)*10),'.mat'),'-struct','S'); 
end